Dynamics of Gas-Fluidized Granular Rods 
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We study a quasi-two-dimensional monolayer of granular rods fluidized by a spatially and tem- 
porally homogeneous upflow of air. By tracking the position and orientation of the particles, we 
characterize the dynamics of the system with sufficient resolution to observe ballistic motion at 
the shortest time scales. Particle anisotropy gives rise to dynamical anisotropy and superdiffusivc 
dynamics parallel to the rod's long axis, causing the parallel and perpendicular mean squared dis- 
placements to become diffusive on different timescales. The distributions of free times and free 
paths between collisions deviate from exponential behavior, underscoring the non-thermal character 
of the particle motion. The dynamics show evidence of rotational-translational coupling similar to 
that of an anisotropic Brownian particle. We model rotational-translation coupling in the single- 
particle dynamics with a modified Langevin model using non-thermal noise sources. This suggests 
a phenomenological approach to thinking about collections of self-propelling particles in terms of 
enhanced memory effects. 

PACS numbers: 



o 



I 

-a 
c 

o 
o 



in 



oo 
o 



X 



I. INTRODUCTION 

Flocking birds, schooling fish, and swarming bacteria 
are examples of collections of self-propelled particles - 
particles that take in energy from their environment and 
then dissipate it by moving in a preferential direction 
through that environment - that display collective, co- 
herent behavior over a huge range of length scales. A 
wealth of theoretical work, ranging from minimal rule- 
based models 0, to coarse-grained hydrodynamic the- 
ories that employ symmetry considerations [E 0, H[ , has 
been conducted in order to elucidate those behaviors 
that are universal to collections of self-propelled particles. 
The earliest physical model [l| predicted a phase transi- 
tion from a disordered to a true long-range-ordered state 
in two dimensions, in direct contrast to thermal systems 
where such a transition is explicitly forbidden. For apolar 
particles with nematic order, the broken-symmetry phase 
is characterized by anomalously large number fluctua- 
tions [E 0] • Other modes of collective behavior predicted 
include swirling and vortex motion and propagating 
waves [IE El El- 

Despite the abundance of models and simulations, lit- 
tle experiment on physical systems has been conducted, 
with what has been done focusing on vertically-vibrated 
granular systems [E EE EE EE EE E3 ■ These systems 
afford an advantage over biological systems: there is 
control over the microscopies and energy input, and a 
steady state can typically be reached provided the sys- 
tem does not age. Vibrated-bed experiments have ob- 
served swirling vortices [E Ell' pattern formation [l4l |. 
and shape-dependent long-range ordering [l5j]. Re- 
cently, giant number fluctuations were observed in a 
vertically-vibrated monolayer of nematic-ordered gran- 
ular rods 1 1 61 ] - Large number fluctuations were also re- 
ported for a collection of vibrated spheres [Hj], fueling 
some debate [l9| as to whether the large density fluctua- 



tions are due to the self-propelled nature of the particles, 
inelastic collisions, or a heaping instability [IH observed 
in vibrated media. Additionally, vertical vibration as 
a method of driving has several inherent disadvantages: 
the driving force is temporally inhomogeneous and the 
particle dynamics may vary depending on the phase of 
the oscillation cycle at the time the particle contacts the 
substrate. Further, for smooth plates, a lateral force is 
only generated by particle overlap out of plane; parti- 
cles 'boosted' to high speeds in this manner contribute 
to power law tails in velocity distributions. 

In the interest of discerning universal characteristics 
of collections of self-propelled particles, we consider an- 
other method of driving. Here, we report on a quasi-two- 
dimensional monolayer of granular rods fluidized by an 
upflow of air. The fluidizing airflow is temporally and 
spatially homogeneous, with airspeed less than free-fall 
and with Reynolds number much greater than 1. Because 
of their extended shape, the particles move by scoot- 
ing when one of their ends has lifted off the substrate. 
This results in self-propelled behavior characterized by 
dynamical anisotropy between translations parallel and 
perpendicular to the particle's long axis. To elucidate 
the effect of self-propulsion, we characterize the dynam- 
ics of a single particle and compare it to the behavior of 
an anisotropic Brownian particle [IE HH, [H|. To that 
end, we model our system using the same Langevin for- 
malism as for a Brownian particle with a single modifi- 
cation: particle self-propulsion is included implicitly via 
non-thermal noise terms. 



II. EXPERIMENTAL DETAILS 

The system we study, depicted in Fig. QTa), consists 
of a monolayer of 435 cylindrical plastic dowel pins - 
length L = 0.95 cm, diameter d = 0.24 cm, and mass 
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FIG. 1: (Color online) A monolayer of gas-fluidized bipolar rods. The diameter of the system is 6.75". a) A frame capture of 
the system from the video data after thresholding. The particles are plastic dowel pins - length 3/8", width 3/32", and mass 
0.055 g - occupying an area fraction of 42%. The call-out shows an actual particle, b) A time-trace of one particle (highlighted 
in the left image). The color code corresponds to instantaneous orientation with respect to the horizontal axis of the image. 



m = 0.055 g - fluidized by an upflow of air. The partic- 
ular area fraction <fi = 42% ensures that the particles are 
uniformly distributed across the system. The particle as- 
pect ratio, L/d = 4, is chosen to prevent any long-range 
ordering. Although an effective harmonic potential ex- 
ists due to interactions with the confining walls [HI, [24[ , 
volume exclusion interactions and the turbulent, chaotic 
mixing of individual particle wakes overwhelms interac- 
tions with the boundaries and eliminates this external 
potential. Excluded volume interactions also serve to 
prevent unidirectional whirling of the particles. 

In-plane motion of the rods is excited by a uniform 
upflow of air at speed U = 370 cm/s. This is lower than 
the terminal free-fall speed, so the rods do not levitate or 
fully lift off the plane. It is also low enough that the rods 
do not "chatter" [25|. Instead, one end of the rod lifts 
slightly off the plane and causes the rod to scoot prefer- 
entially in the direction of the tilt. This is most evident 
when only a few rods are present; here, at finite area 
fraction, the mean free distance between rods is small. 
In addition to scooting, the upflow of air also induces 
random short-time motion. Since the Reynolds number 
based on rod length and air speed is large, ~ 10 3 , there 
is turbulence in the form of irregular wakes shed by the 
rods. The shedding frequency /„ is given by a universal 
value of the Strouhal number: St = f v L/U — 0.18; a new 
wake is generated every time the air flows a distance of 
about 5L [26]. Therefore, the rods experience a corre- 
sponding random kicking force that fluctuates on a time 
scale t v = l/f v w 0.018 s 

The apparatus itself and fluidization method are iden- 



tical to those of Refs. (23l . [24[ . The apparatus is a rectan- 
gular windbox, 1.5 x 1.5 x 4 ft. 3 , positioned upright. A 
circular testing sieve with mesh size 150 /im and diameter 
30.5 cm rests horizontally on top. To reduce alignment 
of the particles with the wall, we place the particles in a 
free-standing cylindrical insert, inner diameter of 17 cm 
and thickness of 0.32 cm, at the center of the larger bed. 
A blower attached to the windbox base provides vertical 
airflow perpendicular to the sieve. The upper and lower 
halves of the windbox are separated by a 1-inch-thick 
foam filter between two perforated metal sheets to elim- 
inate large-scale structures in the airflow. The flow rate 
and its uniformity are confirmed via a hot-wire anemome- 
ter. Raw video data of the fluidized particles is captured 
for 10 minutes at 120 frames per second by a digital cam- 
era mounted above the apparatus. Post-processing of the 
video data is accomplished using Lab VIEW. 

Figure El^b) shows a sample time-trace, 
{x(t),y(t),6(t)}, of a single particle. The color code 
denotes the instantaneous orientation of the particle 
with respect to the horizontal axis. The motion appears 
heterogeneous: we note several long stretches where the 
particle is preferentially moving parallel to its long axis. 
Apart from these, the wandering of the particle is like 
that of an isotropic particle undergoing a random walk. 
To obtain this time-trace, we convert each frame of the 
raw data to binary as it is saved to video from buffer. 
Using Lab VIEW'S "IMAQ Particle Analysis" algorithm, 
we locate the centroid of each particle in the thresholded 
image, imposing an upper bound on the allowable area 
of a single particle, and determine the orientation of 



3 



0.1 K 



* io 3 

10 -5 
10 -6 



parallel 
perpendicular 

■angle 



0.1 0.2 0.3 0.4 0.5 
Kinetic Energy, K [ergs] 



0.6 



FIG. 2: (Color online) Distributions of kinetic energies, Ky 
and Ki, for motion parallel to and perpendicular to the par- 
ticle's long axis and Ke, the rotational kinetic energy. The 
average kinetic energy for each components is {Ku) = 0.021 
ergs, (K±) = 0.012 ergs, and (Kg) = 0.016 ergs. The kurtosis 
excess of the velocity distributions are 0.5 for parallel, 0.4 for 
perpendicular, and 0.6 for angle. A Gaussian distribution has 
a kurtosis excess of 0. 



each particle with respect to a fixed horizontal axis. 
Any particles identified above the area bound consist 
of two or more particles that have collided. In order 
to distinguish the individual particles, a series of image 
processing steps known as erosions, similar to Ref. [25| . 
is carried out. FigureQJa) shows the result of the erosion 
process for one frame in which all particles have been 
separated. 

Next, we link together particle positions and orienta- 
tions from frame to frame by finding the minimum dis- 
placement between two center-of-mass positions in sub- 
sequent frames. To ensure correct matches, we constrain 
this with a maximum allowable displacement and rota- 
tion. The resulting time-traces are then smoothed by a 
running average. We estimate error in the position and 
angle data as <T rms /\/N, where a rms is the rms devia- 
tion of the raw data from the smoothed data and N is 
the number of frames in the smoothing window. This 
yields errors of 18 fim and 3.3 mrads in the position and 
angle data, respectively. Finally, to further minimize any 
wall effects, we include only those segments for which the 
particles are at least three particle diameters away from 
the wall. 



III. SELF-PROPULSION 

As noted in Section [Til the sample time-trace, 
Fig. [JJb), shows long stretches where particle orienta- 
tion is aligned with the direction of motion, indicating 
scooting motion. To quantify this self-propelling behav- 
ior, we obtain the distributions of kinetic energies for 
motion parallel and perpendicular to the rod's long axis 
and for the rod's orientation, shown in Fig. [2j We imme- 
diately see that equipartition of energy does not hold: at 
large energies, the distribution of parallel kinetic energies 
is much greater than the perpendicular distribution. The 



parallel component, with average value (K\\) = l/2m(vjj) 
~ 0.021 ergs, is roughly twice as energetic as the per- 
pendicular component, with (K±) = l/2m(v 2 ) ~ 0.012 
ergs. The rotational kinetic energy splits the difference 
between the two. Its average value is roughly the average 
of these two: (Kg) = 1/21 (lo 2 ) ~ 0.016 ergs, where we 
have used I = (m/12)(3(d/2) 2 +L 2 ). This shows that the 
gas-fluidized rods convert energy provided by the upflow 
of air into in-plane motion, preferentially parallel to the 
long axis. They do so at the expense of motion perpendic- 
ular to the rod's long axis. Such microscopic dynamical 
anisotropy - to which the emergence of collective macro- 
scopic behavior has been attributed [l5j - is a universal 
feature of theories and simulations of self-propelled par- 
ticles, even when the anisotropy is not specified as a part 
of the model a priori. We can rightly consider our gas- 
fluidized rods in the context of self-propelled particles. 

We also obtain the average self-propelling velocity, 
V\\ ~ 0.87 cm/s, from (K«). From this and the lab 
frame diffusion coefficients D x and Dg to be discussed 
in Section [Vj we calculate a = v 2 /(D x Dg) = 27, a 
dimensionless parameter that determines the extent to 
which self-propulsion dominates over stochastic fluctua- 
tions [l~2j. For values a > 4, self-propulsion dominates. 
Thus, we are clearly in a regime where self-propulsion 
effects will be readily observable. Furthermore, we can 
quantify the spatial and temporal extent of the 'scooting' 
behavior by calculating an alignment order parameter 
m(t) — cos[6(t) — 4>(t)], the cosine of the angle between 
the particle's instantaneous orientation, 9(t), and the di- 
rection of its instantaneous velocity, <f>(t) . This quantity 
is 1 if a particle moves in the same direction that it is 
pointing and if it moves perpendicular to its orienta- 



tion. Here, 
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autocorrelation (m(t) ■ m(t + r)), we extract a correla- 
tion time of 6.4 s. We quantify the spatial extent of 
the stretches from the value of the mean square paral- 
lel displacement at the correlation time. This gives a 
displacement of 1.6 cm, roughly one particle length per 
correlation time. 

Self-propulsion is a strictly non-thermal phenomenon. 
The shape of the energy distributions in Fig. [2] further 
indicates the non-thermal character of our system. All of 
the distributions deviate sharply from an exponential at 
small energies. This sharp maximum has been observed 
for a single rod bouncing on a vibrated surface ■ This 
is due to the strong correlation mentioned earlier between 
out-of-plane and in-plane motion. Specifically, when the 
particle begins chattering out-of-plane, its in-plane mo- 
tion is significantly reduced, and vice versa. The non- 
exponential form indicates that the velocity distributions 
are non-Gaussian. We confirm this by calculating a non- 
Gaussian parameter - the kurtosis excess, (v 4 ) / (v 2 ) 2 — 3, 
which equals for Gaussian distributions - of the velocity 
distributions. With kurtosis values of 0.6 for the angle, 
0.5 for the parallel, and 0.4 for the perpendicular, none 
of the velocity distributions are Gaussian. 
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FIG. 3: (Color online) Probability distribution of a) times 
between collisions and b) displacements between collisions. 
The mean free time, r c ~ 0.063 s, and mean free path, r c ~ 
0.042 cm, are labeled. 



IV. COLLISIONS 

The anisotropy of particle shape gives rise to self- 
propelled behavior and also alters the excluded volume 
interaction between particles. As a first step toward 
understanding these steric interactions, we compile bi- 
nary collision statistics. We generate distributions of free 
times and free paths between collisions, shown in Fig- 
ure [3j using a graphical approach. We create a particle 
template that is the size of a grain had it not under- 
gone any erosions or thresholding during data postpro- 
cessing. Then, we reconstruct each frame of the video by 
overlaying the template at the center-of-mass of each ac- 
tual particle and rotating it by the particle's orientation. 
By detecting particle overlap, we determine whether a 
given particle has collided with any neighbors in a given 
frame. We then extract the time increment, At, and dis- 
placement, Ar, between collisions for a given particle and 
compile them as probability distributions of free times, 
Fig. Ela), and free paths, Fig. [3Jb). 

Although the tails of the distributions in Fig. [3] be- 
have exponentially, there is substantial deviation at short 
times and paths. As such, we obtain well-defined mean 
values by taking the average of all the data rather than 
from an exponential fit. This yields a mean free time 
between collisions, r c ~ 0.063 s, and a mean free path, 
r c ~ 0.042 cm, marked as the dashed vertical lines in 
the figure. The ratio, r c /r c = 0.67 cm/s, is near the rms 
velocity, 0.76 cm/s. On subsequent dynamics plots, we 
mark the collision time with an orange dashed vertical 
line. 

For thermal, isotropic systems, collisions are indepen- 
dent events; we expect exponential distributions for free 
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FIG. 4: (Color online) Dynamics of rod orientation: a) mean 
square angular displacement MSDe, b) angular diffusion co- 
efficient De = MSDg/(2r), and c) angular velocity autocor- 
relation function. The line with slope 1 in a) corresponds to 
late-time diffusive motion; the line with slope 2 shows short- 
time ballistic motion. The horizontal line in b) is the rota- 
tional diffusion coefficient, Dg= 0.22 s" 1 . The dashed vertical 
line (orange) marks the collision time, r c ; the dash-dot line 
(green) marks the directional memory time, t$ = (2Dg)~ 1 . 
The solid curves are a Langevin model using a non-thermal 
noise source, given by Eqs. (J7J and 1(8$. 



paths and free times. The deviation of the distribu- 
tions in Fig. [3] from exponential behavior at short times 
and paths shape further indicates non-thermal behavior. 
Distributions similar to ours were observed for spheri- 
cal grains vibrated on an inclined plane (28[. There, the 
sharp maximum was attributed to inelastic clustering. 
Here, it may be attributable to local alignment - which 
may make a subsequent collision dependent on prior col- 
lisions - or to dynamical anisotropy. That is, a particle 
may be more likely to collide at short times as it propels 
through its environment. 



V. LAB FRAME DYNAMICS 

Despite being driven far from equilibrium, one or two 
gas-fluidized spheres act as Brownian particles in a har- 
monic trap [231 ] - However, in Section [ill! we showed that 
a rod moves more energetically parallel to its long axis 
and this thermal analogy subsequently fails. Here, we 
further ask how does self-propulsion alter the dynamics 
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of a gas-fluidized rod? We begin by calculating single- 
particle dynamics in the lab frame where particle motion 
is characterized by both angular and translational dis- 
placements. The lab frame and y, correspond to 
the horizontal and vertical axes of the raw images; an- 
gular orientation is measured counterclockwise with re- 
spect to +x. Recall from Section [II] that we only analyze 
data for which a particle is within 3 particle diameters 
from the wall, effectively breaking a single time-trace into 
many shorter time-traces. A typical particle moves across 
the system in less than a minute. Therefore, dynamical 
quantities, such as MSD#, are truncated for delay times 
greater than 40 s. 

FigureHJa) shows the rotational mean square displace- 
ment, MSD e (r) = ([9(t + t) - 8(t)} 2 ), as a function of 
delay time r. We observe ballistic behavior (ex r 2 ) at the 
shortest time scales and a crossover to diffusive behavior 
(oc r) at long times. 

Rotational diffusion, characterized by the rotational 
diffusion coefficient Dg(r) — MSDe/(2r) shown in 
Fig. (Hb), sets a 'directional memory' timescale Tg = 
(2Dg)~ 1 . At times greater than Tg, a particle will have 
forgotten its initial direction and all directions become 
equal. The long-time value of the angular diffusion co- 
efficient, Dg = 0.22 s _1 , is obtained from the plateau 
in Fig. 01b), giving the value Tg ~ 2.27 s, shown as the 
vertical (green) dashed line in Fig. |4j 

Chiral particles subject to external forces, such as 
those arising from the vibration of a substrate [23 ] or from 
air flowing past them as in the current experiment, spin 
in a preferred direction determined by the sign of their 
chirality. Small manufacturing defects impart chirality 
with a random sign and magnitude to some of the rods 
in this experiment. Some of the long-time rise in Dg(r) is 
due to the spinning of some particles throughout the en- 
tire data set. If these particles are removed, an increase 
in slope is still detected, indicating heterogeneous dy- 
namics for non-whirling particles. Just as the time-trace 
in Fig. []Jb) shows regions where particle scooting is in- 
termittent with thermal-like wiggling, plots of 9(t) show 
regions of fluctuating motion intermittent with rapid ro- 
tation. We stress that, if we exclude whirlers, the long- 
time value of Dg, and thus the timescale Tg, does not 
change significantly. 

The angular velocity autocorrelation, Wg[r) — (ui(t + 
t) ■ Lu(t)), plotted in Fig. BJc), shows a two-step decay, 
characterized by a large, positive rebound near r c due to 
collisions, followed by small oscillations at longer times 
before noise dominates. The behavior of Wg(r) at t c al- 
lows us to deduce the effect that collisions have on the 
particle. Typically, a collision results in a negative re- 
bound: the particle recoils in a direction opposite to its 
incident direction. The large, positive rebound in Wg(r) 
suggests that the nature of our collisions is to re-align 
the particle as to its initial direction. 

We compute the translational dynamics in the lab 
frame, shown in Fig.0 There is no discernable difference 
between data along x and y; only data along x is plotted. 
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FIG. 5: (Color online) Dynamics of rod position in the lab 
frame: a) mean square displacement MSD^, b) diffusion co- 
efficient D x = MSD x /(2r), and c) velocity autocorrelation 
function. The horizontal line in b) shows the long-time 
value of the diffusion coefficient, D x — 0.128 cm-s -1 . The 
dashed vertical line (orange) marks the collision time, r c ; 
the dash-dot line (green) marks the directional memory time, 
Tg — (2Dg)~ 1 . The solid curves are a Langevin model using 
a non-thermal noise source, given by the first term in Eq. [13] 



Figure [5fa) shows the lab frame mean square displace- 
ment, MSD^t) = {[x(t + t) - x{t)f). The behavior is 
ballistic at short times and becomes diffusive at longer 
times. This is confirmed in Fig. 02b): the lab frame dif- 
fusion coefficient, D x {t) = MSDj;/(2t), approaches its 
long-time value, D x = 0.128 cm/s, at about 10 s. The lab 
frame velocity autocorrelation function, W x (r) , shown in 
Fig-EJc), has slower-than-exponential decay with a small 
wiggle at t c . Long-time statistics are poorer and the re- 
bound at Tg may or may not be real. 

Equipartition of energy holds roughly between the lab 
frame and the angle. From W^O) and Wg(Q), we find 
that (l/2)m{vl) - 0.017 ergs, and (1/2)1 (u 2 ) - 0.016 
ergs. Remarkably, despite the self-propelled nature of 
the particles, averaging over all possible orientations re- 
sults in bulk behavior that is more or less thermal in 
appearance. This is consistent with the results in Fig. [2[ 
the distribution of kinetic energies in the lab frame is 
nearly identical to that of the angle and the average of 
the parallel and perpendicular distributions. 



6 



VI. COUPLING ROTATION AND 
TRANSLATION 

In this section, we ask how do rotation and translation 
couple for a self-propelling particle? To explicitly visu- 
alize rotational-translational coupling, we calculate the 
dynamics in a 'fixed-angle' lab frame with axes x and y. 
We construct this frame by rotating the coordinates of 
an entire time trace by the initial particle orientation so 
that x and y are, respectively, parallel and perpendicular 
to the initial direction of the long axis of the particle. 
The axes then remain fixed in time. This is equivalent to 
setting the initial orientation of all particles to da = 0. 

For comparison purposes, it is instructive to review 
how rotation and translation couple for a Brownian par- 
ticle [2(J HH, If a particle is not allowed to rotate, 
translational motion is characterized by anisotropic dif- 
fusion - with two diffusion coefficients, Du and D± - 
for displacements parallel and perpendicular to the par- 
ticle's long axis. The two components, it should be 
noted, become diffusive on the same timescale. If the 
particle is allowed to rotate, this anisotropic diffusion 
regime will cross over to isotropic diffusion characterized 
by a single diffusion coefficient, D x — (1/2)(-Dii + D±). 
The crossover timescale is the same 'directional memory' 
timescale discussed earlier in Section |V] Tg = (2Dg)^ 1 . 
Thus, at times longer than Tg, the fixed-angle lab frame 
axes will become random and the dynamics along x and 
y will become equivalent to the conventional lab frame, 
x and y. 

The fixed-angle lab frame mean square displacements 
MSDj(t) and MSD s (r), Fig. g^a), are both ballistic at 
short times. This short-time behavior confirms the dis- 
cussion in Section IIIII that equipartition of energy does 
not hold, with MSD £ (r) twice as large as MSD^ (r). After 
Tg, the fixed- angle MSDs become diffusive and eventually 
converge, showing that the coordinate axes are random- 
ized and the particle has forgotten its initial direction. 
This indicates coupling between rotation and translation 
in our system. However, examining Fig.[6jb), we do not 
observe a fully-developed anisotropic diffusion regime. 
Instead, following the initial anisotropic ballistic regime, 
D x indicates diffusion at approximately 6 s but D y does 
not become diffusive until about 20 s. Isotropic diffu- 
sion occurs near 20 s, when the two diffusion coefficients 
converge. 

We see that coupling for self-propelled particles is dif- 
ferent than Brownian in two significant ways. First, the 
intermediate-time dynamics show that the two compo- 
nents become diffusive on different timescales. Recall- 
ing the ballistic-to-diffusive timescale for the conven- 
tional orientation-averaged lab frame, x becomes diffu- 
sive sooner and y later. This is most likely because the 
parallel component is more energetic than the perpen- 
dicular component; thus, it takes longer for the y com- 
ponent to 'catch up' and become equal to x. Secondly, 
isotropic diffusion occurs at a timescale an order of mag- 
nitude larger than Tg. Thus, we think of self-propulsion 
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FIG. 6: (Color online) Dynamics of rod position in the fixed- 
angle lab frame: a) mean square displacement, b) diffusion 
coefficient Df = MSD//(2r), and c) velocity autocorrelation 
function. At t=0, x is aligned parallel to the particle's long 
axis. The dashed vertical line (orange) marks the collision 
time, t c ; the dash-dot line (green) marks the directional mem- 
ory time, Tg — (2Dq)~ 1 . The solid curves are a Langevin 
model using non-thermal noise sources, given by Eq. Q13[l. 

as a 'memory enhancement' effect. Isotropic diffusion oc- 
curs much later because self-propulsion, in effect, allows 
the particle to remember its initial direction for a longer 
time. 



VII. BODY FRAME DYNAMICS 

While the fixed-angle lab frame is useful for illustrating 
the effects of self-propulsion at short and intermediate 
times, it is useful to consider yet another reference frame 
to capture the effects at longer times. The body frame 
is a set of cLX6S , X and y, which are re-oriented at each 
time step to coincide with the long and short dimensions 
of each rod, respectively. The individual displacements 
are then summed up successively to form a set of time- 
traces, {x(t), y(t)}. Here, 9(t) = at all times and there 
is no coupling in this frame. 

At short times, the body frame MSDs in Fig. [7fa) are 
identical to the fixed-angle lab frame MSDs of Fig.^a). 
At long times, MSDg crosses over to diffusive behavior 
whereas MSDj becomes superdiffusive (ex t 5 / 4 ). Self- 
propulsion gives rise to enhanced diffusion for transla- 
tions along the particle's long axis. As seen in Fig. [7^b), 
the perpendicular diffusion coefficient, D±, has reached 
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FIG. 7: (Color online) Dynamics of rod position in the body 
frame: a) mean square displacement MSD;,, b) diffusion co- 
efficient Db — MSDb/(2r), and c) velocity autocorrelation 
function. The co-rotating body frame axes are always either 
parallel or perpendicular to the particle's long axis. The line 
with slope 5/4 in a) corresponds to superdiffusive motion. 
The dashed vertical line (orange) marks the collision time, t c ; 
the dash-dot line (green) marks the directional memory time, 
Tg = (2De)~ 1 . The solid curves are a Langevin model using 
a non-thermal noise source, given by Eqs. (fT0|) and (fTTj) . 



its long-time value while Du continues to increase. 

The body-frame velocity autocorrelation functions, 
W||(t) and W±(t), shown in Fig. EJc), reveal the long 
memory of the particle motion. The perpendicular com- 
ponent has a slow decay that exhibits the same rebound- 
ing features at r c and Tg as seen in Wg{r) and W x (r). 
The parallel component shows a smooth, slower algebraic 
decay for the entire run indicating remarkably long-lived 
velocity correlations. This is consistent with the long, 
unbroken stretches of scooting motion seen in the time- 
trace image, Fig. QJb). 



VIII. MODEL 

Collections of self-propelled particles are typically 
modeled in one of two ways: establishing a minimal set 
of rules [3, Q or writing hydrodynamic equations, includ- 
ing all terms consistent with relevant symmetries [!, 0, Q . 
Self-propulsion is usually included as a phenomcnological 
parameter. In keeping with such minimal models, we de- 
scribe our system with a Langevin formalism constructed 
along the lines of that for an anisotropic Brownian parti- 



cle. We make a single modification: rather than write a 
self-propelling force, we implicitly include self-propulsion 
with non-thermal noise terms. We begin by constructing 
Langevin equations for the three independent degrees of 
freedom - the angle and the two body frame components. 

The effective harmonic potential of the bed has been 
eliminated and there are no externally imposed forces or 
torques. The remaining forces - arising from interparti- 
cle collisions, interactions with the substrate, and hydro- 
dynamic interactions with wakes - can be considered as 
noise. Thus, our Langevin equations simply relate time 
derivatives of angle and displacement to random-noise 
torques and forces. The equation for the orientation of 
the rod, 6(t), is 



d9 



Ce(t) 



(1) 



where Ce{i) is Gaussian angular noise with zero average 
and variance 



(2) 



From Eq. fl}, Wg(r) is the velocity autocorrelation func- 
tion of Fig. I^c). 

The equations for displacement can be expressed in 
the lab frame Xi a b(t) = {x(t),y(t)} or in the body frame, 
x(t) = {x(t),x(t)}. In the body frame, x(t) and y(t) 
decouple: 



^-=v x {t)=^{t), 



dy 
dt 



(3) 



(4) 



where C\\ (t) and (j_ (t) are Gaussian random noises with 
zero mean and variance 

(C||(*)C||(0) = («||(*)«||(0> = W r ||(r) (5) 



(a(t)a(t')) = (v±(t)v ± (t')) = w ± (T). (6) 

The noises W±(t) and Wn (r) are the body frame velocity 
autocorrelation functions of Fig. \T[c) . 

The Langevin equations derived above apply to both 
equilibrium and non-equilibrium systems. In equilibrium 
systems, the noise fluctuations of Eqs. @, ©, and 
are determined by the fluctuation-dissipation theorem. 
They adopt white-noise forms when rotational and trans- 
lational friction coefficients do not exhibit memory ef- 
fects. Non-equilibrium systems are not restricted by the 
fluctuation-dissipation theorem. Here, we use the ex- 
perimental forms of the angular and translational veloc- 
ity autocorrelation functions to set Wg(r), W±(t), and 
W\\ (r). Since these quantities contain all the information 
about self-propulsion, we are able to implicitly include 
self-propulsion in our model via non-thermal noise. We 
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will then be able to test whether modeling self-propulsion 
as non-thermal noise is sufficient in lieu of specifying an 
actual force acting along the long axis of the particle. 

The form of the angular noise Ce (i) can be determined 
from the two-step decay of Wg (r) , Fig. [4|c) . This sug- 
gests the sum of two exponentials: 



2Ai 



Wg(r) = ZM-e Hr|/ri 

n 



T-2 



) 



(7) 



where a is a real number between and 1. The fit to 
Eq. ([7|), shown as the solid curve in Fig. 0|c), was made 
by constraining Dg to equal its long-time value, 0.22 s -1 , 
rather than using it as a fitting parameter. Although 
unable to capture the sharp rebound caused by collisions, 
Eq. ([7]) provides a good fit to Wg(r). We extract two 
correlation times: ti= 0.018 ± 0.005 s and t 2 = 0.11 
± 0.01 s. The smaller correlation time is identical to 
the vortex shedding timescale r„ calculated in Section HT1 
The value of r 2 is roughly the collision time. 

We also calculate the analytical form of the mean 
square angular displacement: 

((AO) 2 ) = 2D e [aS(t - t', n) + (1 - a)S(t - t', r 2 )] (8) 



where we define 



S(t, T ) = \t\-T(l-e-^). 



(9) 



Using the fit values from Eq. ([7]), we plot Eq. ([8]) as the 
solid curve in Fig. BJa). As a consequence of constrain- 
ing the value of Dg when fitting to Eq. ([7]), the result 
overestimates the value of the MSD at short times. 

The form of the body frame noise, C|i(i) and (±(t), is 
obtained from W||(r) and W±(t), Fig. [JK C ); respectively. 
The slow decay suggests a power law form: 



Wi(r) 



A, 



(10) 



The solid curves in Fig. [JJc) are fits to Eq. (flU]) , yielding 
exponents of a± = 0.99 ± 0.04 and ot\\ = 0.73 ± 0.02. 
These terms implicitly include the two non-thermal ef- 
fects of self-propulsion. First, the magnitudes, A;, con- 
tain information about the energy gap between the two 
components. Secondly, the power law exponents incorpo- 
rate the extended memory effect of self-propulsion. We 
note that a power law decay, r~ d / 2 , is expected for par- 
ticles suspended in a fluid due to diffusive transport of 
momentum through the surrounding fluid 30] . Although 
this is not the case for our study, we highlight it as a 
potential analogy: particles in a viscous medium and 
self-propelling particles both exhibit velocity autocorre- 
lations with extended memory effects. This form is also 
able to capture the superdiffusive behavior of MSDm(t). 
We obtain the following expression for the body frame 
MSDs: 



(l-ai)b, 



l-ati 
■a, *T. 



(11) 



This form shows a crossover from r 2 at short times to 
t 2_Q; at long times. Using the parameters obtained from 



the fit to Eq. (|10p . we plot the functional forms given by 
Eq. (fTl"]) as the solid lines in Fig. [7Ja). 

We now have enough information to construct the cou- 
pled fixed-angle lab frame dynamics. Writing the velocity 
in the fixed-angle lab frame in terms of the body frame 
velocity, 



(12) 



we obtain for the fixed-angle lab frame velocity autocor- 
relation function: 

W f (t,t') = \[W\\(t) + Wt(T)]e-« Ae ^/ 2 5 tJ 

+\[W\\(r) ^(r)]M y (0(O))e-« e ( t ') +e ( t » 2 >/ 2 , (13) 

where 6(0) = and 



cos26»(i) sin 26»(i) 
sin26l(t) -cos26(t) 



(14) 



We use the fitting parameters obtained from our fits 
to Eqs. (J7J and ((TOJ) to generate Eq. (fT3]), plotted as the 
solid curves in Fig. Oc). The model shows that Wx(t) 
and Wy (t) converge at Tg , consistent with a crossover in 
the data. 

We then numerically integrate Wf(t,V) according to 



([Ax(t)} 2 ) = / dh / dhWfi^h) 



(15) 



to obtain the solid curves in Fig. [DJa). The compo- 
nents of the model converge on the same timescale as 
the data. We see that modelling self-propulsion as an 
external noise source with long-lived correlations is suffi- 
cient to reproduce rotational-translational coupling. Our 
model suggests that a phenomenological way to think 
about collections of self-propelling particles is in terms 
of enhanced memory effects rather than explicitly detail- 
ing novel forces and torques in microscopic equations. 

The model also reproduces the angle-averaged lab 
frame dynamics well (Fig. [5]). Averaging our expres- 
sions over all initial angles eliminates the second term 
in Eq. (flB"]): the result is plotted in Fig.[5fc) as the solid 
curve. We integrate according to Eq. (fT5|) to obtain the 
solid curve in Fig. [5£a). The model describes the angle- 
averaged velocity correlations very well out to 1 s, after 
which the data falls off more rapidly. The good agree- 
ment with the fit here confirms the bulk thermal behavior 
of the collection of self-propelling rods. 



IX. CONCLUSION 



((A£ 



™ \2\ 



2A, 



(1-^0(2-^)6; 



■[(a, + hr) 



We have investigated the dynamics of gas-fluidized 
rods. Particle shape anisotropy leads to dynamical 
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anisotropy, characterized by preferential motion paral- 
lel to the particle's long axis. Ours is a model system - 
with the advantage of a temporally and spatially homo- 
geneous driving method - to further investigate universal 
phenomena predicted for collections of self-propelled par- 
ticles. 

In this report, we compared the coupling of rotation 
and translation couple for a self-propelled particle to 
that of an anisotropic Brownian particle. A modified 
Langevin formalism implicitly specifying self-propulsion 
via non-thermal noise describes the dynamics data well, 
capturing rotational-translational coupling at the correct 
timescale. Despite the energy gap between the parallel 
and perpendicular components, the model was able to 
reproduce the loss of directional memory at long times. 
Furthermore, despite the non-thermal behavior of in- 
dividual particles, the bulk angle-averaged behavior is 
nearly thermal. 

Future work will continue to explore phase space in the 



interest of observing collective behavior and spontaneous 
symmetry breaking for denser collections of both bipolar 
and polar self-propelled particles. We are interested in 
whether we can induce collective macroscopic behavior 
by manipulating the boundaries of the system as well. 
We are currently working on characterizing compression 
waves that propagate through denser collections of gas- 
fluidized rods. We hope that, through comparison with 
theoretical models and recent vibrated-bed experiments, 
our system will further shed light on universal behavior 
of collections of self-propelled particles. 



Acknowledgments 

We thank Adam Abate for helpful discussion and for 
programming assistance. This work was supported by 
the NSF through grant DMR-0704147. 



[1] 


T. Vicsek, A. Czirok, E. Ben- Jacob, I. Cohen, and 


[16] 


V. Narayan, S. Ramaswamy, and N. Menon, Science 317 




0. Shochet, Phys. Rev. Lett. 75, 1226 (1995). 




(2007). 


[2] 


G. Gregoire and H. Chate, Phys. Rev. Lett. 92, 025702 


[17] 


A. Kudrolli, G. Lumay, D. Volfson, and L. S. Tsimring, 




(2004). 




Phys. Rev. Lett. (2008). 


[3] 


J. Toner and Y. Tu, Phys. Rev. Lett. 75, 4326 (1995). 


[18] 


I. Aranson, A. Snezhko, J. Olafsen, and J. Urbach, Sci- 


[4] 


J. Toner and Y. Tu, Phys. Rev. E 58, 4828 (1998). 




ence 320, 612c (2008). 


[5] 


J. Toner, Y. Tu, and S. Ramaswamy, Ann. Phys. 318, 


[19] 


V. Narayan, S. Ramaswamy, and N. Menon, Science 320, 




170 (2005). 




612d (2008). 


[6] 


S. Ramaswamy, R. A. Simha, and J. Toner, Europhys. 


[20] 


F. Perrin, J. Phys. Radium V, 497 (1934). 




Lett. 62, 196 (2003). 


[21] 


F. Perrin, J. Phys. Radium VII, 1 (1936). 


[7] 


H. Chate, F. Ginelli, and R. Montagne, Phys. Rev. Lett. 


[22] 


Y. Han, A. M. Alsayed, M. Nobili, J. Zhang, T. C. Luben- 




96, 180602 (2006). 




sky, and A. G. Yodh, Science 314, 626 (2006). 


[8] 


I. Aranson and L. Tsimring, Phys. Rev. E 67, 021305 


[23] 


R. P. Ojha, P. A. Lemieux, P. K. Dixon, A. J. Liu, and 




(2003). 




D. J. Durian, Nature 427, 521 (2004). 


[9] 


I. S. Aranson, D. Volfson, and L. S. Tsimring, Phys. Rev. 


[24] 


R. P. Ojha, A. R. Abate, and D. J. Durian, Phys. Rev. 




E 75, 051301 (2007). 




E 71, 016313 (2005). 


[10] 


Y. Tu, J. Toner, and M. Ulm, Phys. Rev. Lett. 80, 4819 


[25] 


A. R. Abate and D. J. Durian, Phys. Rev. E 72, 031305 


(1998). 


(2005). 


[11] 


R. Simha and S. Ramaswamy, Phys. Rev. Lett. 89, 


[26] 


S. Rajagopalan and R. Antonia, Exp. Fluids 38, 393 




058101 (2002). 




(2005). 


[12] 


A. Baskaran and M. C. Marchetti, Phys. Rev. E 77, 


[27] 


H. Wright, M. Swift, and P. King, Phys. Rev. E 74, 




011920 (2008). 




061309 (2006). 


[13] 


D. Blair, T. Neicu, and A. Kudrolli, Phys. Rev. E 67, 


[28] 


D. Blair and A. Kudrolli, Phys. Rev. E 67, 041301 (2003). 


031303 (2003). 


[29] 


J. Tsai, F. Ye, J. Rodriguez, J. Gollub, and T. Lubensky, 


[14] 


J. Galanis, D. Harries, D. Sackett, W. Losert, and R. Nos- 




Phys. Rev. Lett. 94 (2003). 




sal, Phys. Rev. Lett. 96, 028002 (2006). 


[30] 


B. F. Alder and T. E. Wainwright, Phys. Rev. Lett. 18, 


[15] 


V. Narayan, N. Menon, and S. Ramaswamy, J. Stat. 




988 (1967). 




Mech. (2006). 







